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1. INTRODUCTION AND SUMMARY 


The purpose of this research is to develop and to apply new nonlinear system methodologies to the 
stability analysis and adaptive control of high angle-of-attack (a) aircraft such as the FI 8. The present 
progress report reviews the research of this project over the first year (actually 18 months with a no-cost 
extension). 

Considerable progress is documented on nonlinear adaptive control and associated model 
development, identification, and simulation. Also, it appears that previously derived results for bilinear 
system (BLS) stability [1 ,2] as well as describing functions [3] can be adapted to the Ostroff PIF controller 
and the thrust vectoring component with dead zone. The latter will receive more emphasis in conjunction 
with the PIF and other controls studied here. 

The analysis has considered linear and nonlinear, longitudinal, high-a aircraft dynamics with varying 
degrees of approximation dependent on the purpose as summarized in Table 1 . In all cases, angle of attack 
(a) or pitch rate (q) has been controlled primarily by a horizontal stabilizer (S h ). In most cases studied, a 
li nge r adaptive controller provides sufficient stability. However, it has been demonstrated by simulation of 
a simplified nonlinear model [4] that certain large rapid maneuvers were not readily stabilized by the 
investigated linear adaptive control but were by means of a nonlinear time-series based adaptive control. 
More details of the research competed by this period thus far are reported in Sections 2 and 3 below and 
in previous semiannual reports. With regards to nonlinear simulation programming, it is shown that with C 
language it is possible to improve computation speed by two orders of magnitude over the previously used 
MATLAB. 
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Table 1. Aircraft Models 


Type 

Purpose 

Remarks/Limitations 

1. Linear perturbations at 

a = 5®, 15°, 35®, 60® 

Local control, check of nonlinear 
system, application of well devel- 
oped linear control methodologies 

Local stability 

Only valid for small maneuvers 
Special case of types 2-5 

2. Gain scheduled (non- 

linear function of a) 
from 1 

Gain-scheduled adaptive control 
based on well developed meth- 
odologies 

Simplified description of complex 
system 

Approximate stability 

May have stability problems 
with small number of reference 
states and/or large fast ma- 
neuvers 

3. Volterra series 

a) at reference states 

b) general case 

Nonlinear adaptive control via 
cross-correlation and/or £ priori 
dynamic structure 

stability approximation 

Simplified dynamic description of 
complex system 

Non-orthogonal series approxi- 
mation 

Sufficiency of 2 or 3 kernels 

Large computation time for 
adaptation 

4. Bilinear system 

a) continuous 

b) BARMA 

Nonlinear adaptive control via 
model reference identification 
(NLMRAC) 

Large computation time 

Bilinearizing controllers may be 
more practical than linearizing 

5. Polynomial time series 

Stability approximation 
Simplified dynamic description 

ones 

Polynomial approximation may 
be more accurate but more 
time consuming than linear or 
bilinear approximation 

6. Neural network 

Potential application to adaptive 
control 

Probably less accurate than 4 
or 5 for a given data set but 
accuracy may be more robust 
outside the available data set 

7. Nonlinear ordinary 

differential model 

Accurate approximation to fast 
large maneuvers for “final® design 
and simulation 

Stability 

Neglects flexible modes and 
other complications 
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2. MODELS AND SIMULATION 


2.1 Basic Nonlinear Dynamics 

The dynamic equations of motion are established by a nonlinear six-degree-of-freedom aerodynamic 
model. In general, the aerodynamic-force components referred to as the center of gravity (CG) is denoted 
as (X, Y, Z). The aerodynamic angular moment vector about the CG is given by (Fl_ FM, FN). The thrust 
vector T is represented in body coordinates as (T x , Ty, T z ). Then the four equations with respect to body 
axes become 

u = rv - qo>sin(0) + — + — (2.1 ) 

m m 


Y T v 

v = pw - ru + gcos(0)sin(4>) + — + — * 

m m 


( 2 . 2 ) 


w = qu - pv + gcos(0)cos(4>) + — 

m 



( 2 . 3 ) 


The moment equations with respect to CG and body axes: 

p = C 41 pq + C^qr + C43FN + C^FL + — (F z »Ty - Py,T x ) + — (Py»T z ~ ^z«^"y) &A) 

'ZZ XX 


q = C s1 pr + Cgj (r 2 - p 2 ) + FM + 



( 2 . 5 ) 


r - C 81 pq + C^qr + C^FL + C^FN 


+ ^83 + ^40 (^xa^y ~ ^y~^x) 


( 2 . 6 ) 


'xx 


The Euler equations: 


0 = qcos(<|>) - rsin(4>) 


( 2 . 7 ) 
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(2.8) 


<j> = p + qtan(0)sin(<}>) + rtan(0)cos(<|>) 
where the vector (P y , P UI P 7 ) denotes the position vector from the center of mass (CG) to the aerodynamic 
center (AC) and the vector (P xe> P y ,, P ze ) denotes position vector from the center of mass to the engine 
thrust center. The constants in the moment equations (1 .4-1 .6) are functions of the moment of inertia 
quantities (1^, lyy, l n , and l xz ) as follows: 


c„ - g D - 'D 

(2.9) 

= C^XZ (\ a + l w - lyy)/lxX^ZZ 

( 2 . 10 ) 

C 42 = C* (l^ (lyy - 1^) - lja)/lxX^ZZ 

( 2 . 11 ) 

^43 = ^40*xz/*xx 

( 2 . 12 ) 

^51 = Ozz *xx)/'yy 

(2.13) 

^52 = W'yy 

(2.14) 

^61 = ^40 (*xx Oxx “ ^yy) + ^xz)/lxJzz 

(2.15) 

^62 = ^40^xz Oyy ^zz ^xx)/^xx^zz 

(2.16) 

^63 = ^40^xz/^zz 

(2-17) 


The quantities X, Y, Z, FL, FM, and FN depend on the aerodynamic coefficients C D , C y , C L , C t , C m , C n as 
follows: 


O 

II 

-Ol 

& 

0 

(2.18) 

L = qsC L 

(2.19) 
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X = -Dcos(a) + Lsin(a) 


( 2 . 20 ) 


z = 

FL = 


Y = qsC y 

-Dsin(a) - Lcos(a) 
(qsbC, + P V Z - P Z Y) 


( 2 . 21 ) 

( 6 . 22 ) 


(2.23) 


PM - + Px 2 ) 


(2.24) 


FN = (qsbC n + P X Y - P y X) 

^ zz 

q is dynamic pressure, s effective area, and a, b, c moment arms. 


(2.25) 


Angle of Attack, Sideslip, and Total Speed 

With respect to body axes, the angle of attack, a, the sideslip, p, the total speed, V, are defined as 


a 

= tan- 1 te] 

(2.26) 


UJ 


p 

- sin ” ft) 

(2.27) 

V = 

u 2 + V 2 + w 2 

(2.28) 


Mathematical Structure of Aerodynamic Coefficients 

The mathematical structure of the aerodynamic coefficients are based on the wind tunnel test data 
for the high angle of attack vehicle. The aerodynamic coefficients are considered to be functions of the 
following control variables as well as angle of attack, sideslip, Mach number, altitude, roll, pitch, and yaw 
rates: aileron deflection, rudder deflection, and stabilator deflection. The effects of leading edge flap, trailing 
edge flap, speed brake, landing gear, etc., are not considered. 
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Drag Coefficient: 


Cq - C 0 (a,M,h,6 h ) 


(2.29) 


Lift Coefficient: 


C L = + 


■— pLq (ot,M,h) q + C u (o,M,h) a] 


(2.30) 


Pitching Moment: 


C„ - C mJ («,M.h,5 B ) * ^ [C™, (a.M.h) q * C„. («.M,h) a] 


(2.31) 


Side Force Coefficient: 


C y = C y0 (a,B,M.6 af 5 r ) ♦ C^a.M.hJp 


+ i [(^(a.M.hJp + C^a.M.hJr] 


(2.32) 


Rolling Moment Coefficient: 

C, • C„(a,3,M,8,,5 r ) <■ C„(a,M,h)p 


* A [C„<«,M,h)p ♦ C.(«,M,h)r] 


(2.33) 


Yawing Moment Coefficient: 


C„ = C n0 («,p,M I 6„6 r ,5 h ) + C n( j(o,M,h) p 


+ ^ [C^a.M.hJp + C w (a t M,h)r] 


(2.34) 
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Range of States Variables in Aerodynamic Coefficients 


a (angle of attack) 

-10° to 90° 

p (sideslip) 

-20® to 20“ 

M (Mach number) 

0.2 to 2.0 

h (altitude) 

0 to 60,000 ft 


Control Variables and Their Limits 


3 a (aileron deflection) 

-25“ to 25° 

« r (rudder deflection) 

-30“ to 30“ 

3 h (stabilator deflection) 

-24“ to 10.5“ 

3 t (throttle) 

30* to 131® 


Longitudinal Dynamic Equation 

Assume that the motion of the airplane can be analyzed by separating the equations into two groups. 
The X-force, Z-force, and pitching moment equations comprise the longitudinal equations and the Y-force, 
rolling, and yawing moment equations are called the lateral equations. Longitudinal dynamic equations are 
given by (2.1), (2.3), and (2.5). 

Aerodynamic coefficients are functions of angle of attack, total speed, Mach number, etc. 

We can choose state variables as a, V, q, and 0 instead of u, w, q, and 0 and we assume that v = 0, 
p = 0, r = 0, 4» = 0. 

From relationship between angle of attack (a) and air speed (V) 


u = Vcoso 
w = Vsina 
uw - uw 


-10° < a < 90° 
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From above, we define normal acceleration by 


n z = Vsin(0 - a) + Vcos(0 - a) (0 - a) 

= Vsin(0 - a) + Vcos(0 - a) (q - a) 

This normal acceleration term can be very important in state feedback and in calculating controller gains 
in general. 

2.2 Linear Perturbation Simulation 

The linear perturbation equations were derived by means of a Taylor series of the above model at the 
four reference states corresponding to a = 5° , 15° , 35° , 60® . For the FI 8 data the short period eigenvalues 

are given by A. 1 2 = -0.559 ± j 0.337 (damping C = 0.386, period T = 4.7 sec) and for the phugoid X 34 = - 

0.0085 ± j 0.073 ({’ = 0.117, T = 86 sec). Figures 1 through 8 show example angle of attack and pitch 
rate responses to step changes in horizontal stabilator at 5® and at 60® with and without phugoid mode 
component in the simulation. Neglecting the phugoid is similar to assuming constant air speed. Note that 
the short-period responses (and eigenvalues) agree quite well with NASA simulations [4]. 

2.3 Simplified Nonlinear Models 

The question of interest here was to investigate the possibility of developing an accurate but simplified 
nonlinear model that would remain valid in a large range of operating conditions and at the same time would 
be capable of rendering nonlinear phenomena occurring in high angle of attack post stall flight regime. 
Usual practice in aircraft modeling is to characterize its dynamics by providing so-called stability derivatives 
for different operating conditions. Stalford et al. [5] proposed using a Volterra series approach for 
longitudinal as well as lateral aircraft dynamics, claiming that obtained models characterize the systems 
behavior much better than piecewise linear ones. However, they did not try to find a global model (i.e., valid 
for a large range of angle of attack), developing instead four Volterra series submodels obtained from the 
expansion around for equilibria corresponding to different ranges of angle of attack. Although their 
approximate piecewise Volterra model indeed gives results that agree very accurately with the original model 
derived from wind tunnel experiments data, yet this is hardly due to the inclusion of higher-order terms. In 
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fact, the contribution of the latter to the final solution is almost negligible. The piecewise-linear model 
obtained by omitting higher-order terms from the Volterra approximation is seen to give almost the same 
results as the original nonlinear model, as can be observed in Figures 9-10. The comparison in [5] was 
made with a piecewise linear model obtained from linearization around equilibria different than those used 
for Volterra series expansion, although the approximation of plunging force coefficient resulting from such 
a linearization seemed much better than from first terms of Volterra series expansion. This, however, does 
not mean that nonlinear modeling of aircraft dynamics has no advantages over piecewise linear models. 
The conclusion that can be drawn is merely that while piecewise linear models may accurately predict 
complex nonlinear behavior they are very sensitive to the choice of points of linearization and that the best 
piecewise linear fit to the curve of one of the model coefficients does not necessarily have to give the best 
dynamical model. 

The simplified longitudinal aircraft dynamics mode) described in [5] was taken as a basis for 
investigation of nonlinear phenomena that may occur in high angle of attack regime of flight. The model 
is as follows: 


a = q + 9.168c z (a) - 1.8336(5 h + 7°) + 7.361904 
4 = 5.73(a - 1.55) ♦ 2.865 


(2.35) 


where: 

a = angle of attack (deg) 

q = pitch rate (deg/s) 

5 h = elevator control (deg), horizontal stabilator 
c z (a) = plunging force coefficient (Figure 1 1 ) 

The nonlinearity is seen to come from the angle of attack. The state plane portraits of the system for 
different constant values of control are shown in Figures 12-17. Stable equilibria of the system correspond 
to values of control less than 9.49 or greater than 12.24 or to angle of attack less than 14.74 or greater than 
18.87. In the region between those equilibria an unstable and a stable limit cycle occurs. Interesting 
phenomena can be observed for the zone near to the onset of unstability. As seen in Figures 1 5 and 1 6, 
the equilibrium is still stable and at the same time two limit cycles exist, the inner of which is unstable and 


13 











pmofiiil iorcc coemcient 











separates the areas of attraction of the equilibrium and of the stable limit cycle. So for different values of 
control and state variables the behavior of the system is essentially different Some investigation was made 
about possible ways of characterization of such a behavior by means of a global discrete time nonlinear 
model. A proposed form of the model was as follows 

a(k+1) - £ P w a(k)'q(kyu(k)* 

IJ,< (2.36) 

q(k*1) « E PvaWqfkyuM' 

where the summation is over all possible products of powers that give as a result exponent not greater than 
N. Of course, for a specific model only some of the terms will be taken - the choice based on 6 priori 
knowledge about the systems nonlinearities and the significance of a given coefficient's contribution to the 
approximation. 

In order to identify the model of form (2.36), the experimental data was first collected. The experiment 
consisted of observing the outputs of the system (given by (2.35)) subject to random steps of control. To 
capture such phenomena like limits cycles In the data, the steps were rather long - 40 sec. There were 64 
such steps. The time discretization was chosen to be O.t sec. As a result, the identification data contained 
25,600 points In a state plane for 64 values of control. Then, for a few arbitrarily chosen models of form 
(2.36), the parameters were found by minimization of the quadratic criterion 
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(2.37) 


min (£ (y(k) - mf 

P \ k 

where y(k) stands for a or and y(k) is obtained from (2.36). Among the models tried, the most accurate 
approximation of (2.35) was given by the following one: 

a(k+1) = p 1a o(k) + p 2 „a 2 (k) + p^a^k) + 

P4*qO<) + P 5 « c l( k ) a ( k ) + Pe.PW® 2 ^) + Pt.PW® 3 ^) + 
p 8 «u(k) + p 9 .u(k)a(k) + p 10 u(k)a 2 (k) + p 11lt u(k)o 3 (k) + p 12 « 

(2.38) 

q(k+1) = p 1q a(k) + p^a^k) + p^a 3 ^) + 

P4qP( k ) + P5qP( k )«( k ) + P8qP( k )“ 2 ( k ) + P7qP( k ) a3 ( K ) + 

PsqU(k) + P 9q u(k)a(k) + p 10q u(k)o 2 (k) + p 11q u(k)a 3 (k) + p 12q 

with the following values of parameters 


Pi. 

= 

8.4320*1 O' 1 

Piq 

= 

-5.3094*10‘ 1 

p 2 « 

= 

6.7979*10'’ 

P 2 q 

= 

-1.1410*10’ 3 

Pa* 


-1.2527*1 O' 4 

p 3q 

3 

4.3451 ‘lO -6 

p 4<* 

= 

9.6900*1 O' 2 

p 4q 

a 

9. 7427*1 O' 1 

P&* 

= 

5.81 42*1 0" 4 

p 5q 

= 

-4.7215*10' 5 

P&. 

= 

-5.4326*1 0‘ 5 

p 6q 

■ 

-1 .8024*10'® 

P7 a 

= 

1.3799*1 O' 8 

P7q 

= 

6.3993*10‘ 7 

Pa. 

= 

-2.2902*1 O' 1 

p 8q 

= 

-8.2199*10'’ 

Pa* 

= 

2.9968*1 O' 2 

P9q 

= 

2.3537*1 O' 3 

PlOa 

= 

-1.21 58*1 O' 4 

PlOq 

= 

-5.7795*1 O' 5 

Pn<* 

= 

4.241 0*1 O' 7 

Pllq 

= 

-2.6468*1 O' 7 

Pl2a 

= 

-3.9140*10'’ 

Pl^ 

= 

3.4865*10'’ 


Although some of these values seem negligible, it should be noticed that with the values of a going to 20 
the terms multiplied by the coefficient in question become of order 10 5 which makes their contributions 
significant enough. The identified model (2.38) was tested by calculating its responses for the same initial 
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conditions and values of control as used to get state space portraits on Figures 12-17. The comparison 
between the original data obtained by simulation of (2.35) and the identified model behavior is shown on 
Figures 18-23. It can be seen that limit cycles are accurately rendered by the model, as well as the stable 
zone behavior, although large discrepancies occur with the control values close to the stable/unstable zones 
border. These inaccuracies may be due to insufficient identification data and/or improper choice of 
nonlinear terms in (2.36). This can be helped by testing the hypotheses about the significance of every 
particular coefficient based on residuals with and without it. This procedure could be performed once for 
a given aircraft and the values of resulting set of parameters could then be updated based on on-line 
identification during the flight itself. 

The model (2.38) was then used for developing a nonlinear controller for (2.35). 
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Figure 19 . Original vs. Discrete Time Identified Model 
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Figure 23. Original vs. Discrete Time Identified Model 
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3. ADAPTIVE CONTROL APPROACHES 
3.1 Nonlinear Volterra-Based Control 

As the Votterra series representation of the nonlinear plant dynamics is a natural generalization of the 
linear system characterization by impulse response or by transfer function, it seems natural to incorporate 
the concepts of the Volterra series model into the nonlinear control. Although a lot has been published on 
mathematical theory of Volterra series, including existence, realization, and (to a much smaller extent) 
identification, very little has been done to apply Volterra series to control. In some papers (e.g., [6]), Volterra 
series serve only as a conceptual starting point from which a switch follows to discrete time nonlinear time 
series. Among a few that attempt to build Volterra series controllers, the majority of them deals with discrete 
time systems. The controllers proposed are mainly predictive ones where the control is obtained by solving 
the Nth-order polynomial equation [7,8]. The model is given in the form of discrete time kernels which 
require a tremendous amount of data. The kernels are, of course, truncated at certain time value which 
results in characteristic jump in the control and the output step response after the time corresponding to 
the truncation. 

Continuous time controllers based on Volterra series were systematically developed in [9] with 
formulae for the controller’s kernels given those of the plant and of the desired feedback system. In 
particular, the problem of so-called exact feedback linearization was solved here. However, those formulae 
are of limited practical value because of the properties of Volterra series under feedback. The problem is 
that even finite (e.g., second-order) Volterra series of the open loop results in infinite Volterra series of the 
closed loop. This makes it necessary for the controller to include theoretically an infinite number of 
compensating terms even for a quadratic system. The same problem for the discrete time systems was 
treated in [10]. Instead of time kernels, they used multidimensional Z transforms and they arrived at the set 
of formulae equivalent to those in [9], However, they provided also a very elegant transformation of the 
exact linearization problem solution which results in a controller requiring only as many Volterra terms as 
there are in the controlled plant. The control system obtained in [10] is shown below. 
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H 1h -JI m are the multidimensional discrete transfer functions of the controlled plant and D, is a linear 
controller designed only for the linear part of the plant, i.e., 

0,(2) • G„(z)/(H,(z)(1 - G^z)) (3 ') 

where G r#f is the desired transfer function from U c (z) to Y(z). Thus, with exact knowledge of the nonlinear 
part of the plant’s dynamics, the whole design reduces to the choice of the linear controller using any of well 
known methods. The inner feedback loop will compensate for all nonllnearitles. One of the drawbacks of 
this method Is the necessity of inverting the linear part of the plant However, if only the linear and nonlinear 
parts exhibit the same time lag this operation does not represent a major problem. In the configuration of 
controller shown above finding a control amounts then to solving Nth-order polynomial equation in u. Of 
course, the problem arises whether and how many solutions exist to this equation. In the case of multiple 
solutions a rule of thumb would be to choose the one within the operating range of control values. In the 
lack of solutions (e.g., due to inaccurate modeling), a possible remedy could be, for example, to take the 
real part of a complex root. (In a few simulations performed with very inaccurate models this method 
worked surprisingly well.) 

One very attractive feature of this controller is that its structure makes it possible to utilize it not only 
with models represented in the form of Volterra series, but in fact with any model with easily divided linear 
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and nonlinear parts of the dynamic equations. Let us take, for example, a nonlinear time-series model of 
the form 


y(k) = L(y(k- 1 ),.~.y(k- M),u(k- 1 )„„.u(k-M)) + 


A/(y(k-1),~.,y(k-M),u(k-1),...u(k-M) 


(3.2) 


where L is the linear operator and N Is the nonlinear part consisting of higher-order terms (i.e., with all first- 
oider partial derivatives over u and y vanishing at zero). Now the structure of the controller will be as 
follows: 



This diagram Is equivalent to the following algorithm for the calculation of the control value at the moment 
k. 

a) calculate the output of the linear controller u L (k) 

b) calculate the predicted value of the output at the moment k 

m - /.(y(k-1)....y(k-M).u(k-1)^,u(k-M)) 
/V(y(k-1),...,y(k-M),u(k-1), M ..u(k-M)) 

c) solve the equation for x(k) 


W( 9 (k).y(k- 1 ),~,y(k-M+ 1 ),u L (k)- x (k).u 0 c- 1 )— , u (k-M+1)) - 
- L(jt(k)^f(k-1),...,^(k-M+1),^(k),y(k-1)^.,y(k-M+1)) 


(3.4) 





d) calculate the control 


u(k) = u L (k) - x(k) 

It is worth noting that in such a discrete time realization there is one more feedback loop interaction than 
shown on the diagram. The models in the inner feedback loop do not use the previous values of output 
estimates 9(k-1 ), y(k-2), ... (as would be necessary in the continuous time case) but the real measured values 
of output. So it becomes clear that the above algorithm becomes a sort of prediction controller which tries 
to estimate the effects of the previous controller which tries to estimate the effects of the previous controls 
knowing the previous values of outputs and then to adjust the current value of control so that the nonlinear 
part of predicted output is canceled. 

Nonlinear Control of the Longitudinal Aircraft 

The discrete time nonlinear control algorithm presented above was used for angle of attack 
stabilization and control of the nonlinear longitudinal aircraft model described in Section 2.3. For the 
purpose of controller design, the model (2.38) was used with the parameter set derived from an off-line 
identification process. The linear controller was designed for the linear model with controlled output 
chosen to be 

« L ( k +1) ■ Pi.«lW + P4«q L ( k ) + Pa« u ( k ) 
q L ( k+1 ) - Plq®t_( k ) + P*,q L ( k ) + PeqUW 

The design was performed to obtain the closed loop behavior of the form 

G(z) = 0.05/(z 2 - 1.6z + 0.65) 

In order not to cancel the zero of the plant, the observer polynomial (z-0.7) was also introduced. The 
algorithm for the control value u(k) is as follows. First the estimate of the output at moment k is calculated. 
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«(k) = P 1a «i + P 2 ,«i + P3««i + 

Pu* + P 5 .<1i a i + Ps.^®? + Pr.Qi°i + 

Pa« u i + P9« U 1«1 + Pio. u i«i + Pii. u i®i + Pi 2 « 

4(k) ■ Plq«1 '+ P2q«? + p3q«1 

P 4 q®ll ♦ PsqQl ®1 + P8qPl“l + P7q<ll®1 + 

Psq^l + P 9q U 1 a 1 + Pi 0 q U 1®1 + PllqUl®? + Pl2q 

where a 1 = a(k-1), q, = q(k-l), u 1 = u(k-i). 

The linear portion of this estimate will be 

«lW = Pi«“i + p4«^i + P«« u i 

$L<k) = Plq«1 + P4q<ll + Psq U 1 

and the nonlinear portion 

i N (k) = a(k) - « L (k) 

4„<k) - 40<) - ^u(k) 

Then the nonlinear portion of the output estimate in the moment k+1 given control u(k) is equal: 

a N (k+1) = a (k+1) - a L (k+1) 

= Pi.(«-®o + P2.° 3 + P4.M0 + Ps-di + p 6 .q« 2 + p7.q« 3 + 

P9« U + PlO« U “ 2 * Pl1« U “ 3 + Pl2« 

where a = a(k), a L = a L (k), q = 3(k), q L = Q L (k), u = u(k). 

The control value comes from the equation 

« N (k+1) = Pi.®#) + P 4 AM + Pa«(u L (k) - u(k)) 
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which finally yields a solution 


(3-6) 

with «(k) and $(k) taken form (3.5). It is seen that if there are no nonlinearities in the model the control 
reduces to a regular linear controller u = u L . 

A number of simulations was run to test the controller performance, especially in the unstable range 
of angle of attack. Figures 24-28 show the response of the system to the step change of the setpoint of 
angle of attack. The resulting trajectories are compared with desired trajectories following form the linear 
controller design. It can be noticed that modeling inaccuracies do not achieve prefect model following but, 
nevertheless, the system is successfully stabilized and the transients are very smooth and without significant 
overshoots. By different choice of the reference model it is possible to obtain much faster, but at the same 
time much more "nervous" transients. The elevator control as shown as an example on Figure 28a is also 
relatively smooth and, worth noting, its values doe not at all come out from the range corresponding to the 
terminal equilibria. This cautiousness of the controller is the main reason for rather slow regulation process. 
It also can be noticed that some kind of linearizing the closed loop system was indeed accomplished 
because the shape of trajectories is very similar regardless of the zone in which the regulation takes place. 
The reaction of the system for an input disturbance in the form of an impulse of magnitude -1 8 additive to 
the control (i.e., sudden displacement of elevator) is depicted in Figures 29-31. The performance is not 
astonishingly good but still the task of stabilization and disturbance rejection is successfully fulfilled. Of 
course, purely linear controller constant on the whole operation range is not able to stabilize and control 
the plant as can be seen in Figures 32-33. The linear model used for its design was obtained by 
identification from the same data as in the case of model (2.38). 

Conclusions 

The conclusion that comes from the above simulation experiments is that it is possible to model 
nonlinear aircraft dynamics in the form of nonlinear discrete time model containing a limited number of 
power nonlinearities. The proposed nonlinear controller structure was shown to give quite satisfactory 


P Sa u L (k) - (P 2 «« 2 + Ps.« 3 + P5,«q + P»«3« 2 * P7«$« 3 * P 12 .) 
(Ps. + P9.« + Pl0.« 2 + Pll.® 3 ) 
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Figure 24. Step Response with Nonlinear Controller vs. Nominal Response 



Figure 25. Step Response with Nonlinear Controller vs. Nominal Response 



Figure 26. Step Response with Nonlinear Controller vs. Nominal Response 
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results even with the model that was far from perfect. This is, however, a purely empirical result and 
theoretical investigation of such properties of the discussed controller like stability sensitivity of modeling 
errors would be In order. Furthermore, a more accurate discrete time nonlinear presentation of nonlinear 
aircraft dynamics is necessary. 

All simulations included in this section were obtained using PC-MATLAB and its Runge-Kutta 
integrating routine. The minimization of model square error for purpose of identification was performed 
using NekJer-Mead nonlinear simplex method coded in C and compiled by TURBO-C compiler. 
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3.2 Nonlinear MAC Algorithm 

Model algorithmic control (MAC), described, for example, in [11], consists in general in solving the 
model equation for the value of control necessary to obtain required value of output. Usually this desired 
output trajectory is generated from the setpoint by means of a reference model. In case when this model 
is linear the algorithm in essence becomes a linearizing one. More precisely, the equation to be solved is 

Y,*(k + 1) = YmodO^I) + (y(k) - Ymod(k)) (3 - 7) 

where y ref is a desired reference output and y mod is the prediction of the output based on the model of the 
plant. The correction term (y(k) - y mod (k)) takes into consideration the possible error of the model and, in 
fact, introduces integral action into control. In a situation when y(k) is not yet available at the time when 
the control u(k) is computed, as is often the case due to time delays and/or time needed for solving (3.7), 
the correction term may be taken as (y(k-1) - y m0 d(* < - 1 )). anc * ymod(* < + 1) must be based on the 
measurements from moment k-1 which means that the algorithm becomes two-steps-ahead. 

In of model (2.35) with the controlled output assumed to be the angle of attack the algorithm 

takes the form: 


with 


= amod( k+1 ) + («( k ) - «mod( k )) 


«mod( k+ 1) = P. T <M k ) 


(3.8) 


(3.9) 


<J>(k) = [«, a 2 , a 3 , q, qa, qa 2 , qa 3 , u, Ua, Ua 2 , Ua 3 , l] T (k) < 3 10) 

As the control at the moment k must be already computed at moment k the values of a(k) and q(k) are not 
available for its computation so their estimates must be used instead. The correction term is taken to be 
the prediction error from the moment k-1 and the equation becomes 

«„(k*1) ■ * («(k-1) - «^(k-1)) < 3 "> 
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with 


«mod( k+1 ) = P. T *( k ) 

$(k) = [a, a 2 , a 3 , q, qa, qa 2 , qa 3 , U, Ua, Ua 2 , Ua 3 , l] T (k) 

5(k) = P. T <fr(k-1) + («(k-1) - « mod (k-1)) 

m - Pq 4>(k-1) + (q(k-1) - Qmod(k-I)) 

The controller is assumed to know the values of angle of attack and of pitch rate at the moment k-1 . Then 
it estimates their current values a(k) and q(k) taking into consideration previous prediction errors and based 
on them it calculates the control required to achieve a ref at the moment k+ 1 . The value of control is found 
as: 

U(R) = 5 r - Pi«« ~ P2.& 2 - Ps.« 3 - P4«3 ~ Ps,3« - Pa«3° 2 - P 7 «3« 3 - Pi U (3 12) 

P&» + lp 9 «* Pio«« 2 + Pii«« 3 

where 

a r = «^(k + 1) - («(k-1) - a mo< j(k-1)) (3-13) 

and a = a(k), C\ = q(k) as described above. 

This algorithm was simulated for the plant (2.35) with model (2.38) and its parameter values. The 
results of the simulations are seen in Figures 34-35. The reference trajectory was chosen to be 
1 /(z 2 - 1 ,6z + 0.65). The actual output of the plant is seen to follow the reference very closely, even though 
the region of operation was that of the most severe nonlinearities. The control action is also remarkably 
smooth. It should be pointed out that for all simulations presented here setpoints of a correspond to 
equilibria with some negative pitch rate and in reality would result in some decrease of pitch angle, which 
is not included in the model (2.35). Thus, the conditions simulated are somewhat fictitious from the 
aeronautical point of view. Nevertheless, for the purpose of evaluating the performance of control strategies 
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they may provide useful information about the behavior of the system in nonlinear regime. However, more 
detailed study is required. 


3.3 Adaptive MAC 

The discrete time nonlinear state space model (2.36) describes the behavior of the plant (2.35) quite 
accurately in the entire region of operation. Very often, however, such a global model is rather difficult to 
fit and, consequently, one should look for local approximations, depending on the current operating 
conditions. In such a situation, adaptive control seems to offer an ideal solution. For a given model 
structure the controller identifies its local parameters and appropriately adjusts its action. At the same time 
it can compensate for the changes of the plant “true" global parameters. Linear adaptive control has been 
used for nonlinear systems trying to modify the first-order approximation depending on the operating 
conditions treating the nonlinear system as a time-varying linear one. In general this approach requires that 
the plants parameters change slowly, which excludes the case of rapid maneuvers. Nonlinear adaptive 
control is believed to be a proper solution to this problem. While it may be difficult to find a suitably simple 
global approximate model, local behavior may be still highly nonlinear thus causing the linear control to fail. 

The algorithm discussed in the previous section was made to be adaptive, or self -tuning, by 
incorporating on-line identification of the parameters. A recursive least squares (RLS) algorithm was 
implemented in the form taken from 


p(k) = 


Q(k-2) <Kk-1) 

A.(k-1) + <|>(k-1) T Q(k-2) <|>(k-1) 


(3.14) 


Q(k-l) 


1 

Mk-i) 


Q(k-2) 


Q(k-2) 4>(k-1) <t>(k-1) T Q(k-2) 
A(k-1) + <|>(k- 1 ) T Q(k-2) 4>(k- 1) ^ 


(3.15) 


e(k-1) = y(k) - p T 4»(k-1) 


(3.16) 


where y may denote a or q and p may stand for p a or p q , respectively. The forgetting factor X was 
introduced to enable the algorithm to change the estimates of parameters with the change of operating 
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conditions. To avoid the unlimited growth of covariance matrix Q at the steady state when the input is not 
persistently exciting, the variable forgetting factor policy was implemented: 


m = i 


e 

e(k ) 2 


(3.17) 


where e(k) is the current prediction error, e(k) is the acreage prediction error from last 10 samples and e 
is equal to 0.01 . As an additional precaution the trace of the covariance matrix Q was monitored and Q was 
reset to diagonal matrix whenever the threshold value was exceeded. Starting values of parameters were 
taken to be as in (2.36). 

Figures 36-37 show the simulation results of the above algorithm for the same reference trajectory and 
initial conditions as discussed in the previous section. Figures 38-43 display the simulation results for 
another reference model specified as 1/(z 2 - 1.8 + 0.82). Remarkably exact following of the reference 
trajectory may be observed, although, surprisingly enough, the performance is slightly worse than in the 
nonadaptive case. Most probably, this is due to the fact that prediction error now changes much more 
quickly because of the ongoing identification process. Thus, approximating the term (y(k+ 1) - y mo<;1 (k+ 1)) 
by (y(k-l) - y mod (k-1 )) may worsen the behavior of the system as two values of y mod no longer correspond 
to the same parameter vector. Since the on-line identification process assures (at least in principle) that the 
prediction error should asymptotically converge to zero, it is possible that the correction terms in a(k), q(k) 
and in control equation (3.11) ought to be omitted. This will be soon verified in proper simulation 
experiments. 

The performance of the adaptive nonlinear MAC controller was compared to the linear one, which 
uses the same control strategy but with strictly linear model being identified and used for the calculation of 
the control action. The simulation results are shown in Figures 47-51. The initial conditions and reference 
trajectories were exactly the same as the ones for the nonlinear case in Figures 36-44. The starting values 
of model parameters were taken from off-line identification over the entire region of interest (similarly to 
those of model (2.35)). Clear difference between the performance of linear and nonlinear controller can be 
seen in Figures 40-41 and 48-49, particularly in control action at the setpoint a = 15®. The linear identifier 
has 
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Figure 39. Nonlinear Adaptive MAC 
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Figure 43. Nonlinear Adaptive MAC 
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obvious difficulties with fitting the parameters of a linear mode) to the behavior of the plant which is clearly 
nonlinear in this region (see Figure 47). As a result, the control starts oscillating for a while. Also, on the 
other plots it is seen that the nonlinear algorithm results in control plots that are more smooth, although they 
still contain one-pulse spikes. To eliminate these spikes weighting of the increments of control was 
introduced into the algorithm. The controller is designed to minimize the one step ahead cost function: 


J = (ymod( k+1 ) - y r ( k i)) 2 + p( u < k ) - uO*- 1 )) 2 

with y mod , y r as before. Minimization of (3.18) with respect to u(k) yields 

= (y r - a)b + p uQ^I) 

b 2 + P 


(3.18) 


(3.19) 


where 

a = Pi,« + P 2 «® 2 + Ps.® 3 + P4.P + PsA a + Pe.^® 2 + Pz.Q® 3 + Pi 2 « 
b = Ps« + Ps«® ♦ Pio.® 2 + Pii«® 3 

Obviously, for p = 0 (3.19) reduces to (3.12) while for p = <*> we have u(k) = u(k-1) = const. Results of 
simulations of this algorithm with p = 0.02 and p = 0.05 are shown in Figures 52-53 and 54-55, respectively. 
The trade-off between the accuracy of tracking and control smoothness may be observed. For p = 0.05 
the control contains no one-pulse spikes and, in fact, the accuracy of reference following deteriorates only 
slightly. 

Conclusions 

Model algorithmic control based on an approximate discrete state space model works very well for 
the plant (2.35) with angle of attack as the output. Its adaptive version displays behavior slightly inferior to 
the nonadaptive case, which may be a result of inclusion of a probably unnecessary corrective term in 
(3.11). A control increment term in the cost function makes it possible to obtain more smooth control 
trajectories while retaining satisfactory performance. Previous research failed to find good global input- 
output nonlinear time-series approximation for the plant (2.35), so a state space model was used. It seems, 
however, that 
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locally such approximations should be possible thus allowing for application of adaptive MAC algorithm 
based only on input-output data without using state measurements. Theoretical research as to robustness 
of the algorithm with respect to state or measurement noise also will be conducted and relative to more 
complex airframe simulations. 

3.4 PIF Control 

The proportional plus integral plus filter (PIF) control, which has been developed by Ostroff at NASA- 
Langley, has shown remarkable success in simulations as demonstrated in [4] and more recent work at 
NASA. Right tests of a more refined PIF controller are planned for an experimental high-alpha aircraft. 

Basically, the design incorporates linear optimal (quadratic performance index) control about multiple 
equilibrium conditions with controller gain scheduling in conjunction with a tracked command model. We 
are planning to study a simplified version of this controller for possible stability limitations from certain rapid, 
high-alpha maneuvers. We hope to develop a nonlinear adaptive PIF controller and compare its 
performance with the present NASA PIF controller. An example of this configuration is shown in Figure 56. 
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Figure 56. Nonlinear PTF Control Structure 



4. CONCLUSIONS AND PLANS 

Research in the first phase of a three-year program suggests that nonlinear adaptive control (such as 
those based on nonlinear time series or Volterra kernels) can provide improved stability over that of simple 
linear adaptive control. To make this conclusion more conclusive, however, future research will involve more 
complex airframe models and simulations. Also, more complex linear and nonlinear adaptive controllers, 
including PIF designs, will be studied for stability and performance limitations. Comparisons of the studied 
nonlinear control algorithms and the NASA-Langley PIF design (Ostroff) are planned. 
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APPENDIX A 


MATLAB and C Program Comparisons 



NONLINEAR STABILITY AND CONTROL STUDY OF 
HIGHLY MANEUVERABLE HIGH-PERFORMANCE AIRCRAFT. 


NASA Project Progress Report, 
by 

ZION HALEVY ** 

E.C.E, Oregon State University, Corvallis Or. 


JULY 16, 1991 


1 . SUMMARY : 


The object of this report is to summarize part of the 
programming work done on nonlinear adaptive control model 
of aircraft operating in highly nonlinear regimes with large 
values of angle of attack. 

Two models (Cho + and Stalford* models) were simulated on 
PC-MATLAB and converted to "C" . Comparison of the simulation 
time indicated that performance using the "C" algorithm was 
increased by two orders of magnitude. Thus enabling detailed 
parametric analysis in short time span. The fundamental 
structure of C programs can be also utilized extended to 
future simulator models. 


* * On sabbatical leave from ADA, Israel. 

+ Refers to the full nonlinear longitudinal airframe model (Section 2.1 
of Annual Report) 

* Refers to approximate nonlinear model of Reference 3. 
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3. MATLAB SIMULATIONS . 


3.1 Cho model [2] 

An original MATLAB program using an Euler integration 
scheme (BLIN1.M ,COF.M -written by Cho 1991 based on Cao 
et al 1990) was analyzed and restructured in order to enable 
general purpose use and faster execution time. 

The changes performed on the original programs are as 
follows: 

a) Removal of the close-loop BARMA controller (Cho 1990) 
from the main program model BLIN1.M and writing of a 
general purpose structure controller sub program. 

b) Modification and improvement of the structure of the 
program, in particular modification of the sub program COF.M 
(calculate seven aerodynamic coefficients) which was called 
at each iteration. 

The final programs named ZBLIN1.M and ZCOF2.M are in 
the attached disk (see Appendix 3) . 

3.2 Stalford model . [3] 

A model taken from Stalford 1989 was written in 
PC-MATLAB with the same general purpose structure noted 
above . This work was done jointly with J. Dory. 

The final programs named ZRUN5.M ,ZIN5.M and 
ZMAIN5.M using ODE45.M (Runge-Kutta 4th and 5th order 
integration function for numerical solution of ordinary 
differential equations) are on the attached disk. 

The results from the MATLAB simulation ( Appendix 2) 
indicate that the numerical model is identical to that 
presented by Stalford 1989 . 
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"C" PROGRAMS . 


4. 


4 . 1 Cho model 


The final MATLAB programs (ZBLIN1.M , ZC0F2.M) was 
converted to the "C" language (ZBL13.C ,ZC0F3.C and 
ZBL13.EXE -executable program). 

The final programs are in the attached disk (see 
Appendix 3) . 


4 . 2 Stalford model. 

The final MATLAB programs (ZRUN5.M , ZIN5 .M, ZMAIN5 .M, 
0DE45.M ) was converted to "C" language (ZRUN5.C, 
ZIN5.C, ZMAIN5.C, ODE45.C and ZRUN5.EXE -executable 
program) . 


* xxx. M / xxx. C extensions denote MATLAB and "C" programs 
files respectively. 
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5. Comparison between MATLAB and "C " programs 

Matlab is an interactive program that aid with fast 
performence scientific and engineering numerical 
calculation. 

Matlab allows you to solve numerical and simulation 
problems in a fraction of the time it would take to write 
in other languages like "C" or Fortran. 

Furthermore, problem solutions are expressed in almost 
exactly the same way as they are written mathematically. 

It is a complete integrated system including graphics, 
programmable macros, an interpreter and analytical 
commands . 

"C" and Turbo C++ has fast and efficient compiler that 
enables transfer application programs to other systems. 

The MATLAB is an interpreter and hasn't a compiler like 
"C" ,so the "C" language is more efficient . 

The MATLAB simulations included in this report were 
obtained using PC-MATLAB version 3.2 . 

5 . 1 Cho model . 

Comparison of performance are summarized in Table 1. 
Appendix 1 consists of examples plots from both MATLAB and 
"C" in the same conditions as follow: 

Euler integration step= 0.05 [sec] (the desired accuracy). 
tfinal= 10 [sec] (final value of t) . 
t0= 0 [sec] (initial value of t) . 
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TABLE 1 : Cho model. 


PROGRAM 

CONTROLLER 
INPUT dh 

FLIGHT TIME 
[SEC] 

SIMULATION TIME 
[SEC] 

1. BLIN1.M , COF.M 
(original prog.) 

VARIABLE 

(Cho) 

10.0 

1068 

2. BLIN1.M, ZCOF2 .M 
(Improve COF.M) 

VARIABLE 

(Cho) 

10.0 

580 

3 . ZBLIN1 . M , COF . M 

(Remove controller 
from main prog, 
original COF.M) 

CONSTANT 

-1.0 

[deg] 

10.0 

555 

4 . ZBLIN1 . M , ZCOF2 . M 

(Remove controller 
from main prog, 
improve COF.M) 

CONSTANT 

-1.0 

[deg] 

10.0 

364 

5. C PROGRAMS 

CONSTANT 

10.0 

6 


ZBL13.C ; ZC0F3 . C 

(Remove controller 
from main prog, 
improve COF.M) 
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- 1.0 

[deg] 



5 . 2 Stalford model. 


Comparison of performance are summarized in Table 2 . 
Appendics 2 consists of examples plots from both 
MATLAB and "C" in the same conditions as follow: 

ODE45 step integration= 0.1 [sec]. 

ODE45 tolerance= l.e-4 [sec] (the desired accuracy), 
tf inal=T_FINAL=13 . 6565 [sec] (final value of t) . 
t0= 0 . 45*T_FINAL [sec] (initial value of t) . 

TABLE 2 : Stalford model. 


PROGRAM 

CONTROLLER 

FLIGHT TIME 

SIMULATION TIM 


INPUT dh 

[SEC] 

[SEC] 

1. MATLAB programs: 

VARIABLE 



ZRUN5.M , ZMAIN5 . M 

STEP 

7.5 

1866 

ZIN5.M . 

(ZIN5.M) 



2. C PROGRAMS: 

VARIABLE 




STEP 

7.5 

6 

ZRUN5.C , ZMAIN5 . C 

(ZIN5.M) 



ZIN5.C 




(and ZRUN5.EXE) 
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6. CONCLUSIONS AND RECOMMENDATIONS FOR FUTURE WORK. 

6 . 1 Conclusions 

Camprison of MATLAB and "C" simulation reveal that 
performance was increased by two orders of magnitude by 
the "C" simulation . 

MATLAB is a "package" for development stage but for 
detailed parametric analysis the "C" performance better. 

6. 2 recommendations . 

a) Using of PC-MATLAB version 3.5 (for 386 computer) for 
first step developing of the simulation. 

We used in PC-MATLAB version 3.2 for our simulation 
but the new version 3.5 (for 386 computer) is faster 
at least four times (it hasn't limits imposed by 16 
bit nature of 80286/8086 ) . 

We used also PC-38 6/ 20 Epson computer with 80287 
numeric coprocessor chip. 

b) Using of MATLAB MEX-files enables to combine "C" and 
MATLAB programs . 

MEX-file produce from compiled "C" linked into .EXE 
files and renames to .MEX extension. 

So it is possible to call your own "C" programs from 
MATLAB as if they were built-in MATLAB function. 

Speed improvement of up to a factor of 25 are possible 
in this way . 

c) Using of "C" (turbo C++ Ver 1.0 ) simulation for 

parametric analysis. 

d) Using at Cho model in ODE45.M Runge-Kutta integration 
function instead of Euler integration scheme to get 
better accuracy of results. 
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User's Guide of Borland 
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9. APPENDICES 


1. Cho model : examples plot ( Fig 1 - Fig 6 ). 

2. Stalford model : examples plot ( Fig 7 - Fig 10 ). 

3. Programs description of the attached disk. 
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APPENDIX 3 


The programs description in the attached disk : 

1 . Cho model T21 . 
a) MATLAB programs . 

BLIN1.M- Original main program model (including the Cho 
close loop controller) [2]. 

COF.M - Original function called by main program (BLIN1.M) 
to calculate seven aerodynamic coefficients [2] 

ZBLIN1 . M-The final main program after the changes 
performed on BLIN1.M (see 3.1). 

ZC0F2.M- The final function after the changes performed on 
COF.M (see 3.1). 

"C 11 programs . 

ZBL13.C- The final MATLAB program ZBLIN1.M converted to "C" 

ZC0F3.C- The final MATLAB program ZC0F2.M converted to "C" . 

ZBL1 3 . PR J - " C " project (link) file for ZBL13.C . 

ZC0F3.H- "C" header file for ZC0F3.C (prototype). 

ZBL13 . EXE-The final executable program for Cho model. 

RES. RES- The output results file (alpa, V, q, teta, dh , t) . 

ZLC1.M -The program used for graphic program (at MATLAB) . 
to load RES. RES. 

ZPCC5.M -The graphic program called by ZLC1.M. 

ZPC51.M -The graphic program called by ZLC1.M. 

2. Stalford model (31. 
a) MATLAB programs . 

ZRUN5.M- The macro running program that give initial values 
to the model function (ZMAIN5.M) and also graphic 
program (zp5.m). 
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ZMAIN5.M-The function describes nonlinear aircraft 

equations in the form to be used by ODE45.M 
(integrating routine) called by ZRUN5.M. 

ZIN5.M -The function describes input controller called by 
ZRUN5.M , the output of the function is the value 
of the control (the elevetor angle) . 

ODE45.M -Integration of a system of ordinary differential 

equations called by ZRUN5.M .The function is using 
4th and 5th order- Runge-Kutta formulae . 

ZP5.M -The graphic program called by ZRUN5.M. 

b) "C" programs . 

ZRUN5.C- The final MATLAB program ZRUN5.M converted to "C" 

ZMAIN5 . C-The final MATLAB program ZMAIN5.M converted to 
"C" . 

ZIN5.C- The final MATLAB program ZIN5.M converted to "C" 

ODE45.C- The final MATLAB program 0DE45.M converted to "C" 

ZRUN5.PRJ-"C" project (link) file for ZRUN5.C . 

ZMAIN5.H- "C" header file for ZMAIN5.C (prototype). 

ZMAIN51.H- "C" header file for ZMAIN5.C (decleration) . 

ZIN5.H- "C" header file for ZIN5.C (prototype). 

ODE45.H- "C" header file for 0DE45.C (prototype). 

UTIL. C-Utility program for ODE45.C . 

UTIL. H- "C" header file for UTIL.C (prototype). 

ZRUN5 . EXE-The final executable program for Cho model. 

RES. RES- The output results file( alfa, V, q, teta, 
delta_h, tn) . 

ZLC.M -The program used for graphic program (at MATLAB) 
to load RES. RES. 

ZP5.M -The graphic program called by ZLC.M. 

ZP51.M -The graphic program called by ZLC.M. 
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APPENDIX B 
Project Publications 



PROJECT PUBLICATIONS 


(Supported Wholly or In Part by NASA Grant) 

1. R.R. Zakrzewski, R.R. Mohler, “On Nonlinear Model Algorithm Controller Design,' Proceedings, IFIP 
Conf. Sys. Modeling and Optimiz., Zurich, 1991 (to appear). 

2. R.R. Mohler, V. Rajkumar, R.R. Zakrzewski, “Nonlinear Time-Series Based Adaptive Control 
Applications,* Proceedings, IEEE Conf. Decision & Control, Brighton, 1991 (to appear). 

3. R.R. Mohler, Nonlinear Systems: Vol. 2 Applications to Bilinear Control, Prentice Hall, Englewood 
Cliffs, NJ, 1991. 

4. R.R. Mohler, V. Rajkumar, R.R. Zakrzewski, "On Discrete Nonlinear Seif-Tuning Control," Proceedings, 
Korean Control Conf., Seoul, 1991 (to appear). 



